home *** CD-ROM | disk | FTP | other *** search
- /*****************************************************************************
- YUVPAK - a program for
- FRACTAL IMAGE COMPRESSION
- Color Version
- image.tga ==> image.ifs
- version 2.0
- *****************************************************************************/
- #include <stdio.h>
- #include <stdlib.h>
- #include <graphics.h>
- #include <math.h>
- #include <time.h>
- #define number_flips 8
- #define levels 2
- #define max_patch 8
- #define max_scale 1.2
- long int usual (unsigned char Image[max_patch][max_patch], int size);
- int mapping;
- int main(int argc, char **argv)
- { /* begin main block */
- FILE *in, *out, *cf;
- char *inf, *outf, rmsstr[13];
- int ysize = 256,xsize = 256, min_patch = 4;
- unsigned char Image[ysize][xsize][3], blur[ysize-1][xsize-1], plotdx, plotdy;
- unsigned char Range[number_flips][max_patch][max_patch], Domain[max_patch][max_patch];
- unsigned char DX[xsize*ysize], DY[xsize*ysize], colormap[256][3], YUV, YYUV, UYUV, VYUV;
- int YUVTEMP,reverse;
- long int Domain_Class[levels][ysize - min_patch + 1][xsize - min_patch + 1][2][2];
- long int Range_Class[levels][64][64][2];
- int i, rx, ry, dx, dy, x, y, besti, bestdx, bestdy, patchsize, xx, yy;
- long int class1, class2, class3, class4, count;
- int inlevel;
- long int offset, best_offset;
- long int sumr, sumd, sumrd, sumr_sq, sumd_sq;
- long int Domain_sums[2][xsize - min_patch + 1][ysize - min_patch + 1];
- float fsumr, fsumd, fsumrd, fsumr_sq, fsumd_sq, fmagica;
- float best_scale, scale, root_mean_sq, best_root_mean_sq, min_variance;
- float mean_sq, root_mean_sq_tolerance, mean_sq_tolerance, best_mean_sq, fpatchsize_sq;
- float temp, variance, UF, VF, bf, gf, rf, uf, vf;
- time_t start_time, finish_time;
- long time_used;
-
- struct header_t { /* "should" be a 12 byte header... we'll see */
- long time; /* 4 bytes for compression time in seconds */
- short rms; /* 2 bytes for 100.*rms value */
- short add1; /* 2 bytes to be added later... room for growth */
- long add2; /* 4 bytes to be added later... room for growth */
- } header[1];
-
- struct ifs_t {
- unsigned char dx;
- unsigned char dy;
- signed char scale : 7;
- short int offset : 7;
- unsigned short int flip : 3;
- unsigned short int size : 1;
- } ifs[1],
- ifs_table[levels][64][64];
-
- if (argc < 4) {
- printf("usage: yuvpak rms infile.ext outfile.ext \n\n");
- printf("YUVPAK Version 2.0, Copyright (C) 1993, WD Young\n");
- printf("YUVPAK comes with ABSOLUTELY NO WARRANTY\n");
- printf("Please see files 'copying.wy' and 'copying' for details\n");
- printf("If these files are missing,\n");
- printf("write: WD Young, P.O. Box 632871, Nacogdoches TX 75963-2871\n");
- return 1;
- }
-
-
- root_mean_sq_tolerance = atof(argv[1]);
- header[0].rms = (short)(100.*root_mean_sq_tolerance);
-
- inf = argv[2]; outf = argv[3];
- if ((in = fopen(inf, "rb")) == NULL) {
- fprintf(stderr, "Cannot open input file.\n");
- return 1;
- }
- if ((out = fopen(outf, "wb")) == NULL) {
- fprintf(stderr, "Cannot open output file.\n");
- return 1;
- }
- fclose(out);
-
- start_time = time(start_time);
- min_variance = mean_sq_tolerance = root_mean_sq_tolerance*root_mean_sq_tolerance;
-
- GrSetMode(GR_default_graphics);
-
- for (y = 0; y < 64; y++)
- GrSetColor(y,4*y,4*y,4*y);
- for (y = 64; y < 256; y++)
- GrSetColor(y,y,0,0);
-
- for (y = 0; y < 18; y++) fgetc(in);
-
- for (y = 0; y < ysize; y++)
- for (x = 0; x < xsize; x++) {
- bf = fgetc(in);
- gf = fgetc(in);
- rf = fgetc(in);
- Image[y][x][0] = (unsigned char)(0.30*rf + 0.59*gf + 0.11*bf);
- uf = 0.62*rf - 0.52*gf - 0.10*bf;
- vf = -0.15*rf - 0.29*gf + 0.44*bf;
- Image[y][x][1] = (unsigned char)(255.*(uf + 158.)/316.);
- Image[y][x][2] = (unsigned char)(255.*(vf + 112.)/224.);
- }
-
- fclose(in);
- for (y = 0; y < ysize; y+=4)
- for (x = 0; x < xsize; x+=4) {
- UF = VF = 0;
- for (yy = 0; yy < 4; yy++)
- for (xx = 0; xx < 4; xx++) {
- UF += Image[y + yy][x + xx][1];
- VF += Image[y + yy][x + xx][2];
- }
- Image[y>>2][x>>2][1] = (unsigned char)(UF/16.);
- Image[y>>2][x>>2][2] = (unsigned char)(VF/16.);
- }
-
- for (YUV = 0; YUV < 3; YUV++) {
- min_patch = 4;
- if (YUV >= 1) {
- xsize = ysize = 64;
- mean_sq_tolerance = 0;
- min_patch = 4;
- }
- for (y = 0; y < ysize; y++)
- for (x = 0; x < xsize; x++) {
- GrPlot(x,y,Image[y][x][YUV]>>2);
- GrPlot(x+300,y,Image[y][x][YUV]>>2);
- }
- for (y = 0; y < ysize - 1; y++)
- for (x = 0; x < xsize - 1; x++)
- blur[y][x] = (Image[y ][x ][YUV]
- + Image[y ][x+1][YUV]
- + Image[y+1][x ][YUV]
- + Image[y+1][x+1][YUV])>>2;
-
-
- sprintf(rmsstr,"%4.1f",sqrt((double)mean_sq_tolerance));
- GrTextXY(20,290, "rms tolerance", 255, 0);
- GrTextXY(130,290,rmsstr,255,0);
-
-
- /*************************************************************************/
- /* Classify Range's */
- /*************************************************************************/
-
- inlevel = 0;
- for (ry = 0; ry < ysize; ry+=8) {
- GrLine(0, ry, ysize - 1, ry, 384);
- for (rx = 0; rx < xsize; rx+=8)
- {
- for (y = 0; y < 8; y++)
- for (x = 0; x < 8; x++)
- Range[0][y][x] = Image [ ry +y ] [ rx +x ][YUV];
- class1 = usual(Range[0], 8);
- Range_Class[inlevel][ry>>3][rx>>3][0] = class1;
- }
- GrLine(0, ry, xsize - 1, ry, 384);
- }
-
-
- /*
- ****************************************************************
- * *
- * Compute Domain_sums array *
- * *
- ****************************************************************
- */
- inlevel = 1;
- patchsize = min_patch;
- for (dy = 0; dy < ysize - 2*patchsize+1; dy++) {
- GrLine(300, dy, 300 + xsize - 1, dy, 384);
- for (dx = 0; dx < xsize - 2*patchsize+1; dx++) {
- Domain_sums[0][dy][dx] = Domain_sums[1][dy][dx] = 0;
- for (y = 0; y < 2*patchsize; y+=2)
- for (x = 0; x < 2*patchsize; x+=2) {
- Domain[y>>1][x>>1] = blur[dy+y][dx+x];
- Domain_sums[0][dy][dx] += Domain[y>>1][x>>1];
- Domain_sums[1][dy][dx] += Domain[y>>1][x>>1]*Domain[y>>1][x>>1];
- }
- variance = (Domain_sums[1][dy][dx]- Domain_sums[0][dy][dx]*Domain_sums[0][dy][dx]/16.)/15.;
- if (variance > min_variance && YUV == 0) {
- class1 = usual(Domain, 4);
- Domain_Class[inlevel][dy][dx][0][0] = class1;
- Domain_Class[inlevel][dy][dx][0][1] = mapping;
-
- for (y = 0; y < 2*patchsize; y+=2)
- for (x = 0; x < 2*patchsize; x+=2)
- Domain[y>>1][x>>1] = 255 - blur[dy+y][dx+x];
- class1 = usual(Domain, 4);
- Domain_Class[inlevel][dy][dx][1][0] = class1;
- Domain_Class[inlevel][dy][dx][1][1] = mapping;
- }
- else Domain_Class[inlevel][dy][dx][0][0] =
- Domain_Class[inlevel][dy][dx][1][0] =-1;
- }
- GrLine(300, dy, 300 + xsize - 1, dy, 384);
- }
-
- /*************************************************************************/
- /* Classify Range's */
- /*************************************************************************/
- for (ry = 0; ry < ysize; ry+=4) {
- GrLine(0, ry, xsize - 1, ry, 384);
- for (rx = 0; rx < xsize; rx+=4)
- {
- for (y = 0; y < 4; y++)
- for (x = 0; x < 4; x++)
- Range[0][y][x] = Image [ ry +y ] [ rx +x ][YUV];
-
- class1 = usual(Range[0], 4);
- Range_Class[inlevel][ry>>2][rx>>2][0] = class1;
- }
- GrLine(0, ry, xsize - 1, ry, 384);
- }
- /*
- ****************************************************************
- * *
- * Compute Domain_sums array *
- * *
- ****************************************************************
- */
- patchsize = max_patch;
- inlevel = 0;
- for (dy = 0; dy < ysize - 2*patchsize+1; dy++) {
- GrLine(300, dy, 300 + xsize - 1, dy, 384);
- for (dx = 0; dx < xsize - 2*patchsize+1; dx++) {
- sumd = sumd_sq = 0;
- for (y = 0; y < 2*patchsize; y+=2*min_patch)
- for (x = 0; x < 2*patchsize; x+=2*min_patch) {
- sumd += Domain_sums[0][dy + y][dx + x];
- sumd_sq += Domain_sums[1][dy + y][dx + x];
- }
- variance = (sumd_sq - sumd*sumd/64.)/63.;
- if (variance > min_variance && YUV == 0) {
- for (y = 0; y < 2*patchsize; y+=2)
- for (x = 0; x < 2*patchsize; x+=2)
- Domain[y>>1][x>>1] = blur[dy+y][dx+x];
-
- class1 = usual(Domain, 8);
- Domain_Class[inlevel][dy][dx][0][0] = class1;
- Domain_Class[inlevel][dy][dx][0][1] = mapping;
-
- for (y = 0; y < 2*patchsize; y+=2)
- for (x = 0; x < 2*patchsize; x+=2)
- Domain[y>>1][x>>1] = 255 - blur[dy+y][dx+x];
-
- class1 = usual(Domain, 8);
- Domain_Class[inlevel][dy][dx][1][0] = class1;
- Domain_Class[inlevel][dy][dx][1][1] = mapping;
- }
- else Domain_Class[inlevel][dy][dx][0][0] =
- Domain_Class[inlevel][dy][dx][1][0] =-1;
- }
- GrLine(300, dy, 300 + xsize - 1, dy, 384);
- }
- /*
- ****************************************************************
- * *
- * Range Loop *
- * *
- ****************************************************************
- */
- inlevel = 0;
- for (patchsize = max_patch; patchsize >= min_patch; patchsize/=2, inlevel++)
- for (ry = 0; ry < ysize - patchsize + 1; ry+=patchsize)
- for (rx = 0; rx < xsize - patchsize + 1; rx+=patchsize) {
-
- GrLine(rx, ry, rx, ry + patchsize, 384);
- GrLine(rx, ry ,rx + patchsize, ry, 384);
- GrLine(rx + patchsize, ry, rx + patchsize, ry + patchsize, 384);
- GrLine(rx, ry + patchsize, rx + patchsize, ry + patchsize, 384);
-
- sumr = sumr_sq = 0;
- for (y = 0; y < patchsize; y++)
- for (x = 0; x < patchsize; x++)
- {
- Range[0][y][x] =Image [ry +y] [rx +x][YUV];
- Range[1][y][x] =Image [ry+patchsize -1 -x] [rx +y][YUV];
- Range[2][y][x] =Image [ry+patchsize -1 -y] [rx+patchsize -1 -x][YUV];
- Range[3][y][x] =Image [ry +x] [rx+patchsize -1 -y][YUV];
-
- Range[4][y][x] =Image [ry+patchsize -1 -y] [rx +x][YUV];
- Range[5][y][x] =Image [ry+patchsize -1 -x] [rx+patchsize -1 -y][YUV];
- Range[6][y][x] =Image [ry +y] [rx+patchsize -1 -x][YUV];
- Range[7][y][x] =Image [ry +x] [rx +y][YUV];
- sumr+=Range[0][y][x]; sumr_sq+=Range[0][y][x]*Range[0][y][x];
- }
- fsumr=sumr; fsumr_sq=sumr_sq;
-
- /*
- ****************************************************************
- * *
- * Domain Loop *
- * *
- ****************************************************************
- */
- if ((patchsize < max_patch)&&
- (ifs_table[inlevel-1][ry/(2*patchsize)][rx/(2*patchsize)].offset != (-500>>3))) goto cleanup;
- best_mean_sq = 10000000000.;
- count = 0;
- for (dy = 0; dy < ysize - 2*patchsize+1; dy++) {
- for (dx = 0; dx < xsize - 2*patchsize+1; dx++)
- if ((reverse = (Domain_Class[inlevel][dy][dx][0][0] == Range_Class[inlevel][ry/patchsize][rx/patchsize][0]))||
- (Domain_Class[inlevel][dy][dx][1][0] == Range_Class[inlevel][ry/patchsize][rx/patchsize][0])||
- (YUV >= 1))
- {
- GrPlot(dx+300+(patchsize>>1),dy+(patchsize>>1),384);
- DX[count] = dx; DY[count] = dy;
- count++;
- reverse = 1 - reverse;
- for (y = 0; y < (2*patchsize); y+=2)
- for (x = 0; x < (2*patchsize); x+=2)
- Domain[y>>1][x>>1] = blur[dy+y ][dx+x ];
-
- sumd = sumd_sq = 0;
- for (y = 0; y < 2*patchsize; y+=2*min_patch)
- for (x = 0; x < 2*patchsize; x+=2*min_patch) {
- sumd += Domain_sums[0][dy + y][dx + x];
- sumd_sq += Domain_sums[1][dy + y][dx + x];
- }
- fsumd=sumd; fsumd_sq=sumd_sq;
- fpatchsize_sq = (float)(patchsize*patchsize);
- fmagica = (float)(sumd_sq - sumd*sumd/fpatchsize_sq);
- for (i = 0; i < number_flips; i++) {
- sumrd = 0;
- for (y = 0; y < patchsize; y++)
- for (x = 0; x < patchsize; x++)
- sumrd += Domain[y][x]*Range[i][y][x];
- fsumrd = sumrd;
- if (fmagica != 0.)
- scale = (fsumrd - fsumd*fsumr/fpatchsize_sq)/fmagica;
- else scale = 0;
- if (scale*scale < max_scale*max_scale) {
- scale = (signed char) 63. * scale / max_scale;
- scale = max_scale * scale / 63.;
- offset = (long int)(fsumr - scale*fsumd)/fpatchsize_sq;
- offset = (offset>>3)<<3;
- mean_sq = (fsumr_sq + scale*(scale*fsumd_sq - 2*fsumrd + 2*offset*fsumd)
- + offset*(offset*fpatchsize_sq - 2.*fsumr)) / fpatchsize_sq;
-
- if (mean_sq < best_mean_sq) {
- besti = i; best_mean_sq = mean_sq; bestdx = dx; bestdy = dy;
- best_scale = scale; best_offset = offset;
- }
- if (mean_sq < mean_sq_tolerance) {
- goto gotbest;
- }
- } /* end of conditional */
- }
-
- } /* end of Domain loop */
- }
-
- gotbest:
- for (i = 0; i < count; i++)
- GrPlot(DX[i]+300+(patchsize>>1),DY[i]+(patchsize>>1),384);
-
- sprintf(rmsstr,"%8.4f",sqrt((double)best_mean_sq));
- GrTextXY(560,20,rmsstr,255,0);
- ifs[0].dx = bestdx;
- ifs[0].dy = bestdy;
- ifs[0].flip = besti;
- ifs[0].scale = 63. * best_scale / max_scale;
- ifs[0].offset = (patchsize == min_patch||mean_sq < mean_sq_tolerance) ? best_offset>>3 : -500>>3;
- ifs_table[inlevel][ry/patchsize][rx/patchsize] = ifs[0];
-
- cleanup:
- GrLine(rx, ry, rx, ry + patchsize, 384);
- GrLine(rx, ry ,rx + patchsize, ry, 384);
- GrLine(rx + patchsize, ry, rx + patchsize, ry + patchsize, 384);
- GrLine(rx, ry + patchsize, rx + patchsize, ry + patchsize, 384);
-
- }
-
- if ((out = fopen(outf, "ab")) == NULL)
- {
- fprintf(stderr, "Cannot open output file.\n");
- return 1;
- }
-
- if (YUV == 2) {
- finish_time = time(finish_time);
- time_used = (long)difftime(finish_time, start_time);
- header[0].time = time_used;
- fwrite(header, sizeof(struct header_t), 1, out);
- }
-
- for (ry = 0; ry < ysize>>3; ry++)
- for (rx = 0; rx < xsize>>3; rx++) {
- if (ifs_table[0][ry][rx].offset == (-500>>3))
- for (y = 2*ry; y < 2*ry + 2; y++)
- for (x = 2*rx; x < 2*rx + 2; x++) {
- ifs[0] = ifs_table[1][y][x];
- ifs[0].size = 1;
- fwrite(ifs, sizeof(struct ifs_t), 1, out);
- }
- else {
- ifs[0] = ifs_table[0][ry][rx];
- ifs[0].size = 0;
- fwrite(ifs, sizeof(struct ifs_t), 1, out);
- }
- }
-
- fclose(out);
- }
- GrSetMode(GR_default_text);
- /* All done. Whew... */
- return 0;
-
- } /* end main block */
- #include "wdyusual.c"
-